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Abstract 

We determine the light quark masses from lattice QCD simulations incorporating the electromag- 
netic interaction of valence quarks, using the splittings of charged and neutral pseudoscalar meson 
masses as inputs. The meson masses are calculated on lattice QCD configurations generated by the 
RBC Collaboration for two flavors of dynamical domain wall fermions, which are combined with QED 
configurations generated via quenched non-compact lattice QED. The electromagnetic part of the 
pion mass splitting is found to be — m^o = 4.12(21) MeV, where only the statistical error is 
quoted, and similarly for the kaon, 1.443(55) MeV. Our results for the light quark masses are mjf^(2 
GeV)=3.02(27)(19) MeV, mf^{2 GeV)=5.49(20)(34) MeV, and mf^{2 GeV)=119.5(56)(74) MeV, 
where the first error is statistical and the second refiects the uncertainty in our non-perturbative 
renormalization procedure. By averaging over ibe to cancel C(e) noise exactly on each combined 
gauge field configuration, we are able to work at physical a = 1/137 and obtain very small statistical 
errors. In our calculation, several sources of systematic error remain, including finite volume, non-zero 
lattice spacing, chiral extrapolation, quenched QED, and quenched strange quark, which may be more 
significant than the errors quoted above. We discuss these systematic errors and how to reduce or 
eliminate them. 

PACS : ll.25.-w, 11.25.Yb, 11.30.Cp 
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I. INTRODUCTION 



Electromagnetic (EM) properties of hadrons offer a rich source of interesting and important 
phenomena. The patterns of the mass sphttings between charged and neutral mesons or the 
mass splittings among the octet or decuplet baryons are sensitive to the isospin breaking from 
different up and down quark masses and the EM interactions. It is also known that the width 
difference of p"^ and and the off-shell hadronic light-by-light scattering amplitude, which can 
not be measured experimentally, play an important role in the Standard Model (SM) prediction 
of the anomalous magnetic moment of muon. 

Recent developments made in both hardware and software have advanced lattice QCD close 
to the goal of realistic QCD calculations, and large scale, high precision unquenched simulations 
are becoming available y, yj . Statistical errors on pseudoscalar meson masses are well under 
control, with typical sizes of one-half of one percent, or less. Recalling that charged- neutral 
meson mass splittings are of 0{a) ~ 0(1%), where a is the fine structure constant, it is 
expected that once EM interactions are successfully included, it will be possible to determine 
the up and down quark masses from first principles by using such splittings as inputs. Thus, 
among others, one can check the simplest solution to the strong CP problem, m„=0. 

In this work, we focus on the determination of the light quark masses using the charged 
neutral splittings of the light pseudoscalar meson masses. Following the pioneering work in 
Is! , we introduce EM fields in a non-compact form and combine them with QCD gauge fields to 
realize the combined QCD + QED theory on the lattice. While these earlier works were done 
with unimproved Wilson quarks in the quenched approximation of QCD, we employ domain 
wall fermions (DWF) 0, u\ on configurations with two flavors of dynamical quarks, generated 



by the RBC Collaboration 
errors. 



which allows us to realize the calculation with less systematic 



Our calculation does not contain either effects of the dynamical strange quark or the EM 
interactions of the dynamical quarks. Removing the former approximation is now not a,n 
obstacle^. A cost-effective method to avoid the latter approximation has been proposed in [9|]. 



^ The RBC and UKQCD collaborations are generating 2+1 flavor domain wall fermion ensembles of gauge 
fields, with even lower masses and larger volumes than those used in the present study. We plan to use these 
in future calculations of the kind reported on here. 
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Since photon fields are massless and not confined, finite size effects may be a significant source of 
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HQ. 



In addition 



error. We examine these effects in the vector-saturation model 
to the above, we neglect disconnected quark-loop diagrams in the pion two-point functions as 
the statistical errors associated with them are difficult to control with current methods. We 
will discuss this point in some detail in the following sections and conclude that the omission of 
these diagrams does not affect the determination of light quark masses and has marginal effect 
on the determination of meson mass splittings. 

We introduce the EM interactions on the lattice nonperturbatively such that the effects 
include all orders of the a expansion. Our analysis with several different values of a ranging 
from the physical value 1/137 to about 10/137 shows that the a dependence of the splitting 
is well described by a linear fit over most of this range. This fact encourages a perturbative 
introduction of the EM interactions as an alternative 151] ■ 



The rest of this paper is organized as follows. We first give a brief overview on the combined 
QCD+QED lattice calculation in Sec. [Tll The formula for the pseudoscalar meson mass in 
the presence of isospin violation is given in Sec. Illli After introducing the simulation details 
in Sec. IIVI the numerical results and discussion of systematic errors are given in Sec. |Vl and 
then concluding remarks are described in Sec. IVII Preliminary results of this work have been 



reported in Ref. 



16| . and our related activity toward the lattice calculation of the anomalous 



magnetic moment of muon has been reported in Ref. 15|. Finally, we note that a similar 



calcu 
Ref. 



ation as the one given here, but in the quenched approximation, has been reported in 
17|. 



II. QCD-KQED CALCULATIONS 

In order to understand systematic errors in our lattice study, we recall some basic mate- 
rial of the QCD+QED system. In Sec. Ill Al we discuss the global symmetry of continuum 
QCD+QED before and after the spontaneous breaking due to QCD dynamics. In Sec. Ill Bt 
how our approximation with which the neutral pion is calculated affects its mass is discussed. 
Sec. Ill CI deals with the axial Ward-Takahashi identity in the QCD+QED system on the lattice 
with domain wall fermions. 
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A. global symmetry in the continuum theory 



Here we identify the global symmetry of the continuum QCD+QED system and the related 
Ward-Takahashi identity. For now we set all quark masses to zero. 

In pure QCD, the global symmetry of the fermion action is Gqcd = SU{3)l x SU{3)r x 
U{l)v, where U{l)v represents the baryon number charge and SU{3)l,r are independent flavor 
rotations on the left and right handed hght quark flelds, respectively. U{1)a does not exist 
due to the axial anomaly. The symmetry which survives after spontaneous chiral symmetry 
breaking due to QCD dynamics is the vectorial part of Gqcd, i-e., Hqcb = SU{3)v x U{l)v- 

The QCD+QED system is obtained by introducing the EM interaction of quarks q — 
{u, d, s)-^ with the ^7(l)em gauge potential Aem^ as 

'S'em — J ^jem-^emij. J (l) 
Je^m = Ql'Qer.q , (2) 

where is the electromagnetic current, and Qem is the 3x3 matrix of the electric charges of 
the quarks. 



Qem 



^ I ^ 

-I 
yo 4y 



(3) 



It is convenient to express Qem in terms of the generators {T"}^^^ ... g of SU{3), which can be 
given by T" = ^ with A" the Cell-Mann matrices and satisfy tr (T'^T^) = as 

Qem ^T' + . (4) 

The EM interaction S'em breaks a part of Gqcd- This can be seen explicitly by look- 
ing at the modified Ward-Takahashi identity for the fiavor-nonsinglet axial- vector current, 

d^Al = zeAmM ? r , Qem] I'l^q - ^ tr (gL^^") FZF^r.^. , (5) 

where F^ui^j^i, = dfj^Aej^^y — diyAe„^^ is the field strength of the electromagnetic field, -Fem^i/ = 
^SixuXpF^m, and we still have the light quark masses set to zero. The second term on the right 
hand side arises from the QED chiral anomaly. The first term on the right hand side vanishes 



for a = 3, 6, 7, 8 in the Gell-Mann basis while the second term does so for a = 1, 2, 4, 5, 6, 7, 
and a hnear combination of a = 3 and 8, 



2 ~ 2 2 2 



^ 
1 
yOO-ly 



(6) 



Thus the global symmetry present in this system is 

Gqcd+qed = SU{2yL X Sf/(2)^ X f/(l)e„, X U{l)v , (7) 

where SU{2)'j^ ^ are the subgroups of SU{?>)l^r generated by T^, T^, T'^. QCD dynamics 
breaks Gqcd+qed spontaneously down to the vector-like symmetry 

i^QCD+QED = SU{2)'y X f/(l)en, X U{l)v • (8) 

The exact Nambu-Goldstone (NG) bosons associated with this spontaneous symmetry breaking 
in massless QCD+QED are , and a neutral meson corresponding to the generator T'^, 
i.e., dd — ss. 

Now remember that without the second term on the right hand side of Eq. ([5]) due to the 
QED anomaly, four neutral mesons are massless in the chiral limit as Eq. ([5]) then vanishes for 
a =3, 6, 7, and 8. In other words, the QED anomaly may make the meson corresponding to the 
generator Qem massive. This meson acquires its mass through a diagram consisting of two quark 
triangles connected by two photons, which is C(a^). Since we neglect 0(Q;^)-effects throughout 
this work, we will consider this meson as a NG boson up to this approximation. Neglecting the 
QED anomaly term, the pattern of chiral symmetry breaking in massless QCD+QED reads 

Gqcd+qed = 5f/(2)'^ X U{1)1 X SU{2y^ x f/(l)'^ x U{l)v 

^QCD+QED = 5f/(2)'y X f/(l)e^ X f/(l)v' . (9) 

Here f/(l)'/^^^ are the subgroups of SU{3)l/r generated by Qcm, respectively. The Goldstone 
bosons associated with this pattern are the strange mesons K^, and tts, r/s corresponding 
to T\ T\ 



B. mass of the neutral pion 



In pure QCD, as long as m„ = m^, isospin symmetry remains unbroken, and we need 
consider only the connected (quark) diagram in the neutral pion correlation function since 
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the disconnected ones cancel exactly^. However, once the EM interaction and the u-d mass 
difference are introduced, disconnected diagrams no longer cancel and must be included. Then 
two complications arise; one is the mixing of neutral pseudoscalar mesons and the other is the 
appearance of the two-photon state as the ground state^. Due to the well known computational 
difficulty associated with disconnected diagrams, we neglect their contributions everywhere in 
our current study. Thus we need to know how this approximation affects the pion mass. 
Remember that calculating the connected diagram only is equivalent to calculating the mass 
of the neutral kaon, and the neutral kaon is massless in the chiral limit in the presence of the 
EM interaction as discussed in Sec. Ill A[ Therefore the neutral pion mass computed with only 
connected diagrams does not have a term of 0{a) in the chiral limit. 

For clarity, let tts and vr'' denote the neutral pion in the basis of flavor eigenstates and mass 
eigenstates, respectively. In Sec. Ill Al we discussed that including terms of 0(a), the neutral 
pion corresponding to T^, i.e. vts, is massless when niu = rrid = 0. In general, the mass of 
the physical neutral pion, 171^^0, is different from m-Trg due to the mixing with other neutral 
pseudoscalar mesons. Following the classic current algebra treatment of Tc^ — rj mixing jl8] and 
the TTs — 1]' mixing [19] based on ChPT or a version of that, in pure QCD with m„ 7^ the 
contribution to m^o is estimated to be proportional to (m^ — uidY /ml. Repeating a similar 
treatment for the EM interaction, it turns out that the fact that ^^0, including 0{q) corrections, 
vanishes in the chiral limit does not change. Hence, the mass-squared of 7r° begins at second 
order in the isospin breaking (m„ — irtdY and/or a^. Since we neglect effects of this order 
throughout this work, we can identify "the neutral pion" consisting of only connected diagrams 
with that of 7r°. 

In summary, as long as we discuss the squared pseudoscalar meson mass and neglect 0{a'^) 
and 0{{mu — mrf)^)-effects, "the neutral pion mass" in this work is approximately equal to the 
physical tt*^ mass, and the difference can arise at 0{am). 



^ The interpolating operator for 713 is (127511 — d'y^d) / ^/2^ so in pure QCD the four resulting disconnected 
diagrams in the correlation function (vrsTTa) cancel. Once isospin symmetry is broken, the up and down quark 

loops are no longer equal and hence do not cancel. 
^ Since the two-photon state must carry orbital angular momentum, each photon has non-zero momentum. 

With our lattice setup, the possible minimum energy for this state is about 1.3 GeV. Therefore we could 
have neglected this state even if we had included disconnected diagrams since the possible minimum energy 
is much larger than the NG boson masses. 
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C. axial Ward-Takahashi identity with domain wall fermions 



We consider the QCD+QED system on the lattice with two flavors of domain wall quarks, qi 
and q2, which have mass m^ai,! and mvai,2 and charge Qvai.i and Qva.1,2, respectively. These fields 
are valence quark fields. The flavor non-singlet axial Ward-Takahashi (WT) identity associated 
with the chiral transformation of the valence quark fields is given by 



= ("^val,l + "^val,2) P'^ix) + (mval,l - mval,2) q{x)j5 



2 ' 2 



q{x) 



(10) 



where q = (gi, q2)^ , -^ul^^) is the conserved axial- vector current whose form is the same as that 



in the pure QCD domain wall fermion system 



' em, V 





20| , but the link variables are now replaced with 




'111 



P°-{x) and Jlq{x) denote pseudoscalar densities. While P°-{x) contains physical quark fields 



only, J^g{x) is defined at t 
as in the pure QCD case 



le midpoint of the fifth dimension and written in terms of bulk fields. 



20I ]. X!!-{x) is give by 



W = -^$^[^s(x)(l-7^)f/QCD, 



2 ' 2 



X 



X 



V cii 



+f,(x)(l + 7;.)f/icDM(^-/^) 

'cm/iV-*' " /^/j ~ \y era, fiV^ ~ ) J 



2 ' 2 



(12) 



where the lattice spacing is set to 1, e(s) = 1 for < s < Ls/2 — 1 and —1 for Ls/2 < s < Ls- 
X^{x) vanishes when a = 3, e = or Qvai.i = Qvai,2, so it is the lattice counterpart to the 
first term in Eq. ([5]). The formula analogous to Eq. (ITOl) is written down for three flavor case. 
From the analogy to the pure QCD case in the domain- wall formalism 2l|, it is inferred that 
the EM chiral anomaly arises from J^^ only after sandwiching Eq. ffTOj) between physical states. 
Due to the presence of the EM chiral anomaly, ^^(x) is not conserved for eQvai,! 7^ e(5vai,2 
even in the chiral limit, and hence no NG boson shows up after spontaneous chiral symmetry 
breaking in this two-flavor theory, only pseudo-NG bosons do. However, since we neglect the 
O(a^) contributions in this work, the neutral pion can be regarded as a NG boson. 
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III. FORMULAS FOR PSEUDOSCALAR MESON MASSES AND SPLITTINGS 



A. QCD with rriu = nid 



First we remind the reader of the next-to-leading order (NLO) partially quenched chiral 
perturbation theory (PQChPT) formula for the "kaon-like" pseudoscalar meson mass-squared 
in pure QCD. In this formula, all Nf sea quarks have degenerate mass msea and valence quarks 



have masses set to mvai,i 



m,, 



ml 



ANLo(ml-) 

Mi 



Mi(l + 



-1 



rrid = mud 7^ "^vaL2 

Anlo("^^) 
Ml 



23|. 



(13) 



Nf{Ml - Ml) 
16 



[{Ml - Ml)Ao{Ml) + {-Ml, + MDAoiMi,)] 



(Ls - 2L^)Ml + (L4 - 2L,)Nf M^ 



Mt 



Ml 



Mt 



P 

2 Bq mud, 

Bo {ms + mud), 



2 Ml 



Ml 

TV ' 



AoiM' 



2 Bq msea, 

M^ /M^ 



(14) 
(15) 
(16) 
(17) 
(18) 
(19) 



In the above denotes the renormalization scale of the effective theory, / the decay constant 
in the chiral limit, and Lj the Gasser-Leutwyler low energy constants at that appear in 
the (9(p^)-chiral Lagrangian of QCD. In the above, mp denotes the physical (one-loop) mass 
while Mp denotes the tree level mass. We take Nf = 2 according to our ensemble of gauge 
configurations. 



In the limit m^ mud, the formula for the NLO contribution to the mass- 
like" mesons made of degenerate valence quarks is obtained from the above 



squared of "pion- 



22 



23| 



ANLo(r»^) 
MI 



2 

Wf 

_16 
~P 



TT S 

167r2/2 



+ 



2M2 - M? 



MI 



■MMl) 



(L5 - 2L^)Ml + (L4 - 2L,)Nj Ml 



(20) 
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B. niu 7^ rrid effects 



Incorporation of the isospin breaking effects due to 7^ into the above mass formula is 
straightforward. For kaons, this is achieved by simply replacing rriud in Eq. (|T6l) with m„ for 
or rrifi for K^. As for pions, rriud is replaced in E g. (|T5l) with {mu+md)/2. At this order in ChPT 
a term proportional to {mu — md)'^ also appears 2J], but arises from the disconnected diagrams 
that contribute to the 7r° and f] correlation functions in QCD which we do not compute. So, 
this term is omitted from our fits. 



C. QED corrections to meson masses 



At leading order, the correction to the charged meson mass squared is 



^(Ql-Q2)^ (21) 

where Qi is the charge of valence quark i in units of e. Since the meson masses do not change 
under simultaneous interchange of quark mass and charge, the only possibilities for 0{am)- 
terms are 

(Qi + Q2)'(mi + m2), (22) 
(Qi-Q2)'(mi + m2), (23) 
{Ql-Ql)(m^-m^), (24) 

where is the mass of valence quark i. Combining this and the discussion in Sec. [Tll we 
parameterize the EM contribution to the meson masses as 

Aem(m2.) = 6{Q,-Q,f (25) 
+ ^0 {Qi + Qjf {rrii + m,j) 
+ 5+ {Qi - Qjf {rrii + rrij) 

+ {Ql. - Q]) {m - mj) 

+ f^sca {Qi - Qjf (2 nisea) 
+ 5m... {Qi + Q,f. 

To avoid confusion, Qi always refers to the charge of the quark, not the anti-quark. Factors 
of a and Bq have been absorbed into the definition of the low energy constants (5's). 0{a'^) 
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contributions have been neglected in Eq. fl25l) . The last LEG is purely a lattice artifact induced 
by the finite size of the extra dimension for DWF and vanishes in the limit ^ oo. There 
is also a similar small mixing with the physical term, 6{Qi — Qj)"^, which can be subtracted 
by measuring the EM effects on the residual mass. We will come to this point again when we 
discuss the residual mass and results of our fits in Section |Vl The logarithmic dependence on 
the quark mass has been calculated for the unquenched theory j25| and very recently for the 
partially-quenched theory |26i], to NLO, i.e., including all terms of 0{am) for the case of three 
valence {riyai = 3) and three sea [ugea = 3) quarks. However such an analysis is not available 
yet for Ugea = 2; we thus omit such terms. A full treatment of the splittings to NLO must await 
our upcoming calculation with 2+1 fiavors of sea quarks. 

Note that the masses of the charged mesons 7r+ (ud) and (us) do not vanish for massless 
quarks, and the only terms that contribute to the neutral mesons are the ones with 6o and 

^rrires ■ 

So, 5_, and 6sea. represent linear combinations of low energy constants (LEC's) in the 
0(a;p^)-part of the chiral Lagrangian {e.g., see Ref. [26|). They are determined by fitting the 
numerical data to the form given in Eq. (l25l) . However, since our neutral pion does not contain 
disconnected diagrams, the value of 6o extracted in this work, which we take to be the same 
for all mesons, could be different from the physical one for the 7r° and r] mesons. 



28| approximately 



With the above parameterization, it is easy to see that Dashen's theorem 
holds, and that it is violated at 0{am) which we parameterize as 

A.M = { _ 1, (26) 



EM part 



following Ref. 29|]. Had we kept terms of 0{a'^), they would also violate Dashen's theorem. 



Adding the above contributions from QCD and QED, the squared pseudoscalar meson masses 
are obtained at NLO, 

ml = Ml + Anlo(M2) + Kn.{Ml). (27) 
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IV. SIMULATION DETAILS 



We employ QCD gauge configurations with two fiavors of dynamical domain wall fermions, 
generated by the RBC Collaboration 8|. The lattice size is x T = 16^ x 32 with degenerate 
sea quark masses, amgea = 0.02, 0.03, 0.04 in lattice units, corresponding to ms/2 < nisea < 
where is the physical value of the strange quark mass. In the lattice spacing is determined 
to be 1.691(53) GeV using mp=770 MeV, and hence the physical spatial volume corresponds 
to ~ (1.9 fm)'^. Most of the results presented below were determined on the same ensembles 
used in [8;], about 100 lattices at each sea quark mass, separated by 50 monte-carlo time units 
(see Tab. [T]). The source time of quark propagators is set to tgrc = 0. Additional information 



on the configuration ensemble is given in Tab. [TTl For further details, see 

A non-compact form is adopted for the lattice QED gauge action as in [4|| . We first write the 
action in momentum space, and impose the Coulomb gauge fixing condition plus an additional 
constraint on the vector potential in the time direction (see Appendix A). After diagonalizing 
the kernel, the Boltzmann weight for the photon fields can be written in a Gaussian form. 
We then generate them by the random Gaussian noise method with e=l. The photon field in 
configuration space, Aem,fi{x), is then obtained by inverse Fourier transformation. Details of the 
generation of the QED gauge configurations are given in Appendix El It is worth noting that, 
thanks to the non-compact form of the action and thus the simple generation procedure, there 
is no auto-correlation among the configurations. Since we make the quenched approximation 
for QED, the sea quarks do not couple to photons, and so the fine structure constant a does 
not run with the cutoff scale. 

Exponentiating Aem,/^(a;) and the valence quark's electric charge eQvai, we obtain U{1) hnk 
variables, 

{UemA^))^^""' = e''''5-'^°"^'''('^\ (28) 

Configurations for the QCD+QED theory are then constructed from f/qcd,^(2;) x (t4m,/^(a;))'^'^™', 
which are used in the inversion of the Dirac operator for valence quarks. 

We calculate the two-point correlation functions of the pseudoscalar and vector mesons with 
seven valence quark masses, m^ai = 0.015, 0.020, 0.025, 0.030, 0.035, 0.040, and 0.0446 at 
each sea quark mass and compute all possible degenerate and non-degenerate two-point meson 
correlation functions. The last mass corresponds to the bare strange quark mass as determined 
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in [8|. The electric charges of the valence quarks are Qu = +2/3 for up-type quarks and 
Qd = Qs = for down-type quarks. In a preliminary study we took three values of e, 

which correspond to a = 1/137, (0.6)^/(47r), l^/(47r), to examine the a dependence of the 
meson mass splittings. Having found that the relative statistical errors on our splittings did 
not change with a, we fixed a to its physical value for the main calculations reported here. In 
addition, we found that statistical errors in the splittings are dramatically reduced by averaging 
correlation functions with e = ±1 as then the leading 0{e) noise, which vanishes in ensemble 
averages, cancels exactly on each configuration. The dramatic reduction in error is shown in 
Fig. [1] for a representative case. 

The lattice spacing determined from the charged vector meson mass is found to be a^^ = 
1.659(37) GeV, using the experimental value of mp= 775.8 MeV. This is consistent with the 
value 1.691(53) GeV found in [81. Our statistical error is smaller because in this case we used 
an ensemble of configurations twice the size used in [8| , and for most of the pseudoscalar mass 
splittings reported on here, by taking lattices separated by 20 and 25 monte-carlo units. These 
were blocked together in pairs to account for possible auto-correlations. 

Finally, we remind the reader that we have only calculated connected diagrams for the 
neutral mesons. 



V. NUMERICAL RESULTS 



A. residual quark mass 



We first evaluate the residual quark mass in the usual way, using the midpoint pseudoscalar 
density in Eq. ffTOl) . Since the residual quark mass depends on the EM charge, we determine it 
separately for uu, dd, and ud mesons by averaging the following ratio over a suitable plateau 



for each type of meson 



30, 



3l|, 



(^5".(^)^"(0)) 

(P«(t)P"(0))' 



amUam,) = (29) 



where 0{t) denotes the operator 0{x) summed over the three- volume. 

In terms of a low energy chiral expansion, the ratio of correlation functions in Eq. (1291) is a 
constant at LO and receives corrections at higher order (proportional to rriq, rrig, and so on). 
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We therefore define the residual quark mass m^es to be the value of the ratio in Eq. (1291) when 
niq 0. The chiral limit is then defined as ruq —rrires- A similar strategy works here, but 
the definition of explicit chiral symmetry breaking effects in the low energy effective theory 
must be extended to include 0{a) effects. Because the electric charges break isospin, the most 
general form for the residual mass is 

rriresia) = m,es(mg) + Ci(Qi - Qa)^ + C2{Qi + Qs)^ (30) 

where mresifnq) is the ratio in Eq. (1291) with a = and contains corrections to all orders in m^, 
and Ci and C2 are coefficients of 0{a) which vanish when Lg —* 00. 

In fact, what matters in the following is the difference Arrires = 'n^res(«) — '^res('^g) which 
enters in mass-squared differences of the mesons. Because this difference is calculated from 
highly correlated values of the residual mass, it can be determined very precisely in each case. 
For uu, dd, and ud type mesons, we find 

Amres,uu = 7.11(5) X 10"^ (31) 

Am,,,,,,- = 1.73(1) X 10"5 (32) 
^m^esM = 2.36(4) x 10~', (33) 

where we have simply averaged the differences for all quark mass combinations. Notice that 
the ratio of the uu to dd value is very close to 4 as it must be according to Eq. ( l30l) . However, 
Amres,ud is not equal to the average of the uu and dd values. Using Amj.es,ud and either the uu 
or value, we can determine Ci and C2. We find 

Ci = 1.92 X 10"^ (34) 
C2 = 4.00 X 10"^ (35) 

We ignore the tiny statistical errors in Ci and C2 in what follows below, and note that the 
a = value of rrires determined here, 0.001387(39), is consistent with that found in [8|]. 



B. correlation function 



To extract the pseudoscalar meson masses, we measure the 
function and fit it according to the same procedure detailed in 

( {p'^{t)y p*^(o)). 



crit) 



bllowing two-point correlation 



(36) 
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where 

P'Ht) = Yl *('^' ^) ^5 Qjit, x), p-».^(t) = J2 x) 75 qjit, y). (37) 
X x,y 

We consider the Uyai = 3 case, so q(t,x) = {u(t,x), d(t,x), s(t,x))'^, and we calculate 
the above correlation functions for all possible combinations of degenerate and non-degenerate 
valence quark masses. The pion in pure QCD is also calculated, and denoted by vr^. The pseu- 
doscalar meson masses computed from the pseudoscalar two-point functions are summarized 
in Tabs. IIVI - IVIIII The fit range is the same in each case and for all quark masses, 9-16, and 
matches that used in ^. The value of x^/dof for each fit is less than or about 1. 

We calculate the vector meson masses in a similar manner. The results for degenerate u,d 
quark mass mesons are summarized in Tab. IIX[ The quark mass dependence of the charged 
vector meson is shown in Fig. [2] as an example, where the line denotes the fit to 

mp± = ay + 2by [ruq + mres,ud) , (38) 

using only the rriq = mgea = "^vai data points. Figure [3] shows the quark mass dependence of 
the p^-p° splitting where the mgea = "^vai points are fit to 

mp± — rripo = + 2 d^nig, (39) 

and yield a small but nonzero value in the chiral limit (~ 0.5 MeV). However, it is seen from 
Fig. [3] that if the chiral extrapolation is made in a different way, for example by extrapolating 
r/ivai to the chiral limit first, a different result is obtained. It appears some of the instability 
originates with the choice of fit range for the vector mass. In Fig. H] we show the splitting 
calculated for a uniform tmin = 5 for all sea quark masses. The spread in the splittings has 
cleary decreased. We note that the lower values of tmin do not yield good ^^i some of the 
'nT'sea = 0.03 and 0.04 data, however, so there may be some excited state contamination in the 
splittings. In any case, it is clear that more statistics is needed to resolve the vector splitting, and 
in particular its quark mass dependence. Furthermore we have omitted disconnected diagrams 
in the calculation of rupo, thus the nonzero value for splitting is not conclusive. Note however, 
that at large m^ai where the statistical errors are under better control, for all values of mgea 
the EM splitting corresponds to about 0.5 MeV. Experimentally, the p+ — mass difference 
is consistent with zero |32i]. The results for the parameters appearing in Eqs. (1381) and (l39ll are 
tabulated in Tab. |Xl The vector mass splitting is an interesting case for further study. 
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C. size of 0(Q!^)-corrections 



Before going to the determination of the quark masses, let us discuss the size of 0{a'^)- 
corrections. In Fig. [5] we show the a dependence of — m^o and m^o — fn'^Q ■ Three values 
of the electric charge corresponding to a=l/137, 0.6^/ (Air) and 1^/ {An) are examined. The lines 
denote a linear fit to each set of mass-squared splittings with the constraint that the splitting 
vanishes at a = 0. 

While m^o — ""^^q is well described by the linear fit over the whole range of a, — m^o at 
the largest value is clearly away from the fit line. In both cases the data points for a= 1/137 
and 0.6^/(47r) lie on a line including a = within tiny statistical errors. This indicates 0(q;^) 
effects only become significant for e ~ 1, which is not terribly surprising. 

D. low energy constants and quark masses 

Now we turn to the determination of the quark masses. To this end, we extract the low- 
energy constants (LEC's) defined in Sec. IIII CI by fitting the difference of the square of the 
pseudoscalar mass with the same quantity evaluated at a = 0. The fits are uncorrelated since 
including the full covariance matrix makes them unstable, and it is likely ill-determined besides. 
Fitting all of the data at our disposal (61 masses at each sea quark mass), we find the LEC's 
tabulated in Tab. IXII The mass squared differences and fit for m^ea = 0.02 are shown in 
Fig. [6l Restricting the fit range to mj < 0.03 does not change the values significantly as seen 
in Tab. IXlf though the x^/dof value is reduced. The corresponding a = values of Bq and the 
NLO Gasser-Leutwyler constants are given in Tab. IXIII 

The LEC's in Eq. ( ]25|l are extracted using the mass dependence of the form -|- rrires- 
Then the remaining explicit chiral symmetry breaking effects at order 0{amres) are contained 
in the constants 6 and Sm^es- The value of S^res Tab. IXII should be compared with the 
value of 6*2 determined from Arrires- In particular, one should have Sm^es — 2-BoC2- One sees 
that this combination is roughly two times too large. Since the LEC Smr-es is found by linearly 
extrapolating to the chiral limit, it is possible that the discrepancy arises due to the omission 
of chiral logs. Recall that the logs are 0{mqlogmq) compared to the linear terms and so could 
make a difference. Also note that the value of C2 is roughly an order of magnitude smaller than 
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the physical 6 LEG. Because 6 and Ci have the same charge dependence, it is their sum which 
is extracted from the fit. Thus, we must subtract off the amount 2BqCi from 5 to obtain the 
physical LEG. This amounts to about a 20% reduction from the fitted value of 5. Note that 
the higher order effects in the residual mass also affect the logarithms in Eqs. (fT8|) and (!20!) 
which therefore do not cancel exactly in the meson mass-squared difference with QED switched 
on and off. However, these terms are order (9(m Amres) in the chiral expansion, so we neglect 
them. 

The meson mass-squared splittings show little dependence on the sea quark mass. In Fig. [7] 
we show a representative sample of splittings for each sea quark mass. This is reasonable since 
we have not coupled the sea quarks to the photons, so the difference between ensembles is likely 
due to fluctuations in the gluon fields. In partially quenched chiral perturbation theory for 2+1 



flavors 261] , there exists a term at NLO that couples the valence quark charges to the sea quark 
masses. To account for this possibility, we include the 5 sea term in our fit (see Eq. ( l25l) . and 
notice that only the charged mesons are affected). From Tab. |XT]we see 5 sea is small, consistent 
with zero for the fit range including all masses, and only two standard deviations away from 
zero for the reduced range, which includes only the lightest two sea quark masses. The values 
of the other LEG's are largely unaffected by the 5sea term, except 5, which is reduced in the 
first case and increased in the second. Because the evidence is not strong in our data for such 
a term, we focus on the fits with 5sea fixed to zero in the following. 

Using the results for the LEG's, we determine the quark masses as follows. As inputs, the 
three experimental values of , Tn^j^± , and m|^o are taken. We avoid the 7r° since we have not 
determined its true mass to NLO due to the lack of disconnected diagrams. We then iteratively 
solve the set of equations generated from Eq. ( 1271) for each meson for the three unknowns m„, 
m^, and m^. At this stage, the quark masses are the bare lattice values, including the shift due 
to the residual mass. 

Using the non-perturbatively determined quark mass renormalization constant 1/Zm = 



^5=0.62(4) [33|, we obtain the values of the light quark masses, m^^(2GeV) = Zmirrig + m^es), 
and ratios shown in Tab. IXIIIi The error on Zs reflects a statistical and a systematic error from 
the choice of Aqcd in the range 250-300 MeV, which have been added in quadrature. As the 
bare quark mass range used in the fit to meson mass-squared is reduced from 0.0446 to 0.03, 
the physical quark masses change by about one statistical standard deviation, or less. The only 
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change outside of one standard deviation is for the quark mass ratios which are determined 
very precisely. The ratio of up quark to down quark mass comes out to be close to one-half 
while the strange quark to average up-down quark mass is about 28 to 29, depending on the 
fit range. 

Knowing the physical up and down quark masses, the physical charged-neutral pion mass 
splitting is found from the LEC's in Tab. IXIi Dividing the mass-squared difference by 771^^+ +m.,^o , 
we find 'm-„+ — m^o = 3.89(17) MeV using all data, and 4.12(21) MeV using quark masses less 



than or equal to 0.03, which is somewhat less than the physical value of 4.5936 MeV[34]. The 
part of this difference arising from the up-down quark mass difference was estimated a long 
time ago to be about 0.17(3) MeV jisl and more recently, 0.32 (20) MeV 27!. In addition, there 
are still several systematic errors like non-zero lattice spacing and finite volume effects to be 
addressed in the calculation of the mass splitting, so the level of disagreement is not surprising, 
perhaps even encouraging. It is interesting to note that almost the entire mass difference comes 
from the leading term (~ 98%) since the physical up and down quark masses are so small (the 
omission of disconnected diagrams for the 7r° mass has a very small effect here). Similarly, we 
find the EM part of the kaon mass difference, mx+ — rrixo = 1.443(55) MeV for the resricted 
range (1.441(42) MeV for the whole range). QED effects make the charged kaon heavier, just as 
for the pions. In this case, the leading term also dominates, but the 0{ams) terms contribute 
about 23%. Ultimately, the pure QCD Bq^tjIu — m^) term dominates the physical states and 



makes the neutral kaon heavier by 3.972(27) MeV 32]. 



From Eq. (j26ll the breaking of Dashen's theorem at 0{amq) is found to be Aem = 0.337(40) 
or 0.264(43), using the restricted fit range. This is somewhat smaller than the large Nc estimate 
reported in Aem = 0.85(24). Following 29], Aem would have to be roughly 10 in our 
calculation to render the up quark massless. 



E. systematic error estimate 

We now turn to a discussion of the systematic errors in our calculation. In the calculation of 
the neutral pion correlation function, we ignored the contribution from disconnected diagrams. 



18 



Disconnected diagrams contribute to m^o (and m^) at O^arrig) ^ so this is potentially a sig- 
nificant effect. However, we avoided using the 7r° mass to determine the quark masses, so this 
will not affect those estimates. And as mentioned above, the physical vr^ — mass splitting 
is dominated by the LO contribution since the up and down quark masses are so small. The 
main effect would be to alter the value of the LEG 6^ for the 7r° and rj mesons which we have 
assumed has the same value for all pseudoscalar mesons. 

In any study of the EM interactions, finite volume effects may be significant as the photons 
are massless (and unconfined). In order to get a rough estimate for the size of this effect, we 
examine the finite volume effect to 6 using the vector-saturation model [lO, 11, 12, 1^, 14] as 



an example^. Assuming that the finite size effect purely due to QCD cancels in the difference 
between m.,^+ and m^o, this model estimates the mass difference in a finite volume, A^^em = 
m^i — m^o, to be 

q&T' 

where N is the total number of sites, T' is the first Brillouin zone in the momentum space apart 
from = 0, and 

r - Eg, 

At=0 

q, ^ 2 sin (^) . (41) 
Applying our lattice volume and taking the ratio to that in the infinite volume, we obtain 

A.,em(^ 1.9 fm) • • ^ ' 

Thus we find roughly a +10% increase in 6. We expect a similar size of correction for the 
other (5j's. The shift in 6 could affect the determination of the quark masses. Shifting the 



We have sought a proof that such contributions enter at 0{am^). In fact, one can show the leading contri- 
bution, where only a single photon is exchanged between the quark loops, vanishes. Likewise, any diagram 
where one photon and any even number of gluons are exchanged also vanishes if one considers each quark 
loop to be made from renormalized propagators and vertices only. However, if more general diagrams are 
considered, like the above but including a gluon exchange between two separated quark propagators from the 
same loop, the proof no longer holds. In any case, we expect these contributions to be suppressed, at least in 
the weak coupling limit of QCD. 



More correctly, a parameter representing non- resonance contributions, which turns out to be small is 
set to zero. 
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value of physical 6 — 2BqCi by 10%, we find the quark masses change by less than 1%. Thus 
we conclude that the finite volume effect on the quark masses due to the EM interaction is 
negligible. However, removing this effect will enhance the pseudoscalar mass splitting itself 
significantly. 

We did not take into account the renormalization of the quark masses due to the EM 
interaction. However, the QED part of the renormalization is expected to change the quark 
mass by 0{a) ~ 1%. Since this is well within the statistical and other uncertainties discussed 
already, we ignore this effect. In future calculations it will be a simple matter to include these 
effects directly in the non-perturbative renormalization calculation of Z^- 

The use of the quenched approximation for the QED gauge fields results in a leading error of 
order a a| in correlation functions since the dynamical quarks do not interact via the quenched 



photons. However, as noted in [26| the sea quark charge effects that enter at 0{amsea) can be 
dealt with in two ways. First, the log terms come with known coefficients so can be subtracted 
from the lattice results before fitting for the LEC's. Second, the LEC's which come with sea 
quark dependence cancel out of some mass squared differences. The latter is operable here but 
not the former since we did not fit to the chiral logs as they are not known for the Nf = 2 
partially quenched case. We expect sea quark effects to be small in this present study since 
they were treated explicitly as neutral particles with respect to the EM interaction. 

Finally, we note that changing our bare quark mass fit range from 0.015-0.0446 (less than 
1/2 rris to rris) to 0.015-0.03 has little effect on the LEC's or quark masses. Experience in the 
pure QCD case has shown that this range of quark masses is likely to be beyond the range of 
applicability of chiral perturbation theory Q], and x^/dof does decrease for the restricted range. 
So, being conservative, we take as central values those results determined from the pseudoscalar 
two-point function with a quark mass fit range 0.015 < rrif < 0.03. There is some uncertainty 
introduced by the inclusion of the LEG 6 sea, which we set to zero or found to be small (or zero 
within errors) when left as a free parameter. Since the evidence is uncertain for such a term 
in our data, and the sea quarks are not charged, we stick with the fit with 6sea = to quote 
our results for the quark masses. Thus, from the third row of Tab. IXIIII our final values for the 
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quark masses are 

mf^(2GeV) = 3.02(27)(19) MeV, (43) 

m5^(2GeV) = 5.49(20) (34) MeV, (44) 

mJf(2GeV) = 4.25(23) (26) MeV, (45) 

mf^(2GeV) = 119.5(56)(74) MeV, (46) 

mjmrf = 0.550(31), (47) 

m./m^rf = 28.10(38). (48) 

The first error is statistical, the second from the error on the renormalization constant Z^. 

Note that uncertainties due to the absence of the strange sea quark and finite lattice spacing 
were not considered in the above discussion of errors. We leave these issues for future work 
where they will be addressed explicitly by using 2 + 1 flavor DWF gauge conflgurations at two 
lattice spacings. 



VI. SUMMARY 

In this work, we have determined the electromagnetic splittings of the pseudoscalar meson 
masses by calculating correlation functions in a combined background of QCD+QED gauge 
fields which were, however, generated separately. The gluon configurations came from a re- 
cent two flavor domain wall fermion simulation by the RBC Collaboration while the QED 
configurations were generated in the quenched approximation. 

The highly correlated nature of the calculations allows very small effects from QED to be 



e, or larger than 
1 as well. Here 



observed, even though the naive statistical errors on hadron masses are as lar 
the mass splittings themselves. This was seen in the original calculation 
we have gone a step further by averaging masses computed with ± electric charge on each 
configuration to cancel 0{e) noise on each configuration. This has lead to extremely small 
statistical errors on LEC's and physical ratios like m^j/m^. 

The charged to neutral pion mass splitting was found to be — m^o = 4.12(21) MeV 
(for the restricted quark mass fit range) compared to the experimental value 4.5936 MeV 32 1, 
0.17(3) or 0.32(20) MeV of which is due to the up-down quark mass difference alone jigl . 27|. 



A simple model calculation leads to an estimate of ~ 10% finite volume effect. Similarly, the 
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kaon splitting arising from electromagnetism is 1.443(55) MeV (note that it is positive). We 
emphasize that the calculations reported on here were carried out for the physical value of a, 
and we take these results to be encouraging. In the pioneering work by Duncan, Eichten, and 
Thacker ^ which used the quenched approximation for both QCD and QED and employed 
Wilson fermions at a single coarse lattice spacing, the pion mass splitting was found to be 
4.9(3) MeV. Another more recent quenched calculation finds a value that is somewhat higher 



still 



17|. 



Using the physical pseudoscalar meson masses as inputs, we were able to fix the values of 
the light quark masses, m^, rrid, and rris, including effects of QED and violations of Dashen's 
theorem through 0{amg). These are given in Eqs. (H3|) - (l48l) . 

In this work we have neglected the contributions to the neutral masses arising from discon- 
nected valence quark loop diagrams in the two-point correlation functions because they are 
difficult to compute precisely. This leads to an uncertainty of 0(a^, arriq, {rriu — rn,})'^) in the 
neutral pion mass-squared (see Eq. ([5]) and the discussion in Sec. IIIB|) . The inclusion of these 



diagrams is left for future work. But notice that this type of diagram is necessary only for the 
calculation of the 7r° mass and not the K^, and we did not use 7r° mass to determine the quark 
masses. 

Finally, we also computed the vector meson mass splitting and found it to be quite small, 
~ 0.5 MeV, essentially zero, considering systematic uncertainties in our calculation. This is an 
interesting topic that we will investigate further in future studies. 

The study presented here nicely sets the stage for future work that will focus on the 2-1-1 flavor 
dynamical DWF conflgurations generated by the RBC and UKQCD collaborations 
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36|, and 



those that will soon be generated by the LHPC, RBC, and UKQCD collaborations. Using 
smaller quark masses, larger lattices, and at least two lattice spacings, the accuracy of the 
meson and quark mass splittings will improve signiflcantly. Present and future work also 



includes the electromagnetic splittings of the baryons 34 1 
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APPENDIX A: GENERATION OF NON-COMPACT U{1) GAUGE FIELD CON- 
FIGURATIONS 



Here we describe the way to generate configurations of f/(l) gauge fields on the lattice in the 
non-compact formulation. In the non-compact U{1) lattice gauge theory, the gauge potential 
Aein,n{x) is treated as a basic dynamical variable and put on the mid-point of the link {x, x + 'fl) 
of the hypercubic lattice with topology and with lattice spacing equal to 1. The purely gauge 
action is 

SnCU{1) = XI idi^AmA^) - ^i.'^emA^) j , (Al) 

where 9^ denotes the forward difference operator 

d^f{x)^f{x + il)-f{x). (A2) 

The gauge potential Acni,^{x) is assumed to obey the periodic boundary condition. Then the 
gauge potential is expressed in momentum space as 

An.A^) = ^ E ^'"'^^'-'^^Mp) ■ (A3) 



per 



Here y = JJ^ A^^ with the number of sites along the /x-th direction, and F denotes the first 



Brillouin zone. 



I 27r 



^. = -(§^-l),•••,-l,0,l,•••,^^ (A4) 
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In the decomposition (lA3p . not all of the modes are independent of each other because Ae^^^{x) 
is real-valued. Using the reflection operator R in the first Brillouin zone F 

1^ TT if = TT 

the reality condition is expressed in momentum space as 

(e^^ MRip))) = (e'^ . (A6) 

In terms of these variables, the action flAll) becomes 

1 I ~ ^2 

= J2 \p^A,{p)~p,A^{p) , (A7) 

pgr o<fi<u<3 

where 

p, ^ 2 sin (^) . (A8) 

Let us introduce a single fermion field ipi^) with a unit charge in the system. The following 
discussion applies even when various matter fields with different charges coexist in so far as 
the minimum charge is redefined to be unity and the system is invariant under the gauge 
transformation of the form 

AemA^) ^ = Aen.,^(x) + d^A{x) , (A9) 

^{x) ^ tp'ix) = e*^(^) ^{x) , (AlO) 
^{x) ^ ^'(x) = ^(x) e-'^^^'^ . (All) 

We recall that from Eq. ( ]A9I) the Wilson line constructed from Aem,fi{x) behaves as usual 

g-«^om,,i{a;) _ g«A(x) g-iylem,M(a;) ^-iA(x+p:) ^ (A12) 

The matter fields are coupled to the gauge potential through 6*^°™''*^^) as in the compact lattice 
formulation. We also remark that A(x) parameterizing gauge transformation is not necessarily 
periodic as long as the fields transformed via Eqs. (1A9I) . ( lAlOp and (lAllll continue to satisfy 



respective boundary conditions. The only condition for A{x) required from this consideration 
is 

A(x + N^i2) = A{x) + 27rr^ (r^ G Z) . (A13) 
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The quantization condition for arises from the presence of matter fields with nonzero electric 
charge. In general A(x) satisfying the condition ( 1A13P can be written in the form 

3 

A(x) = 5^27rr,^ + A(°)(x), (A14) 

where A*^'^''(a;) satisfies the periodic boundary condition. A*^''^(a;) thus has Fourier decomposition 

AW(x) = -^5^e^^-A(p). (A15) 
per 



From the reality condition A{R{p)) = yA{p)j , the constant mode A(0) is real. Using Eqs. (]A3p . 
(lAlSp . the gauge transformation ( lAOp for the gauge potential becomes in the momentum space 



Mp) ^ K(P) = Mp) + 2vrr^ ^5p,o + ^Pf.Hp) • (A16) 

A(0) acts only on the matter fields through e'^^^^\ which ranges over a compact space. We can 
thus leave it unfixed in the gauge fixing procedure for the calculation of the expectation values 
of operators. 

The non-compact formulation needs the explicit fixing of the U{1) gauge symmetry. We 
employ the Coulomb gauge fixing condition 

3 

J2d*A,^,j{x)=0, (A17) 
where 9* is the backward difference operator 

(A18) 

In momentum space the condition (lAlTp becomes 



J2p,Ajip) = 0. (A19) 

i=i 

The condition ( 1A17[) is not sufficient to fix all U{1) redundancy; only the redundancy cor- 
responding to the parameters {A(po, p)}p^o, where p = (pi, p2, ps), is eliminated by this 
condition. Accordingly, for p 7^ 0, a component Aj{p) with pj 7^ is determined by the other 
two spatial components 

Mp) = -irJ2pkMp), (A20) 
25 



irrespective of whether pq ^ or po = 0. The residual gauge symmetry is generated 
by {A(po, 0)}po^o (spatially uniform gauge transformation) and r^. They act only on 
{^^=o(PO) 0)}pQ^o and the constant modes ^4^(0) of A^^^^^x) respectively. For po 7^ 0, a 
given Aq{pq, 0) will be converted to a gauge configuration Aq{pq, 0) = on the same orbit via 
A(po, 0) given by 

A(po, 0) = zilo(po, 0). (A21) 
Po 

In this way we can impose a condition 

I^=o(po, 0)=0 (po^O). (A22) 

to eliminate A{pq, 0). The remained ones are G Z. These can be eliminated by imposing 
the following condition on the constant modes 

Vv 

0< A^(0) <27r^. (A23) 

Thus we eliminated harmful gauge redundancy. 

The degrees of freedom (lA23p . the counterpart of Wilson loops in lower dimensions, cannot 
be gauged away. However such degrees of freedom do not play vital roles for the dynamics 
in four-dimensional gauge theory with sufficiently large volume. Thus, we fix A^j(O) to the 
constants in the range (1A23|) as a boundary condition 

MO) = c^. (A24) 



Now we turn to the description of generations of gauge configurations. It is sufficient to 
concentrate on generating configurations in momentum space since the Fourier transformation 
allows to convert them into the ones in coordinate space. 

First we consider a mode with p 7^ 0. Without loss of generality we can then assume that 
Ps 7^ 0. The independent integrated variables are A^{p) (/i = 0, 1, 2) while A^i^p) is given by 
Eq. (]A20p with j = 3. Inserting such A^^p) into the corresponding part of Eq. (]A7|) . we get 



y- 

^ 2e2 



P 



Aoip) 
2 



k=l 



+ 



p' 



ipsY 



PjPkRe (Aj{p) (Ak{p) 



j, k=i 



(A25) 
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where 



fi=0 



(A26) 



The action (]A25P gives the Boltzmann weight in generating configurations for the modes with 
p 7^ 0. However that form is not useful yet for generating configurations because the two 
spatial components are mixed with each other due to the terms 



M 



where a 2 x 2 real symmetric matrix M takes the form 



Mp) 



(A27) 



M = p 



1 + 



[pif P1P2 \ 



iPs? 



P1P2 



(Ps)' 

, iP2] 



(A28) 



V (ps)' {%Y J 

As this is the eigenvalue problem in two dimension, it is possible to resolve this mixing analyti- 
cally. This is an advantage of Coulomb gauge over the covariant Lorentz gauge. The eigenvalues 
of M are 



{Pif + {P2f 



(A29) 



Those two eigenvalues become degenerate if and only ii pi = Q = In this case, we can use 
the basis {Ai{p), A2{p)) to generate a configuration according to the action 



2e 



p ^2 
p 



(Pi = = P2) ■ 



(A30) 



Here hp = 1 for R{p) = p and hp = 2 for R{p) 7^ p, because the action gets doubled due 
to the contribution from the complex conjugate partner in the latter case. In the case that 
iPiY + (^2)^ 7^ 0, M can be diagonalized as 



M = O diag (m_, m+) O . 



Here the orthogonal matrix O is 



O 



r2 ri 
-ri r2 



(A31) 



(A32) 
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with 



(j = 1,2). 



We first generate a configuration for (a^{p), A^{p) ) according to tlie action 

hr, 



(A33) 



2e2 



A-ip) 



(A34) 



A configuration for ( Ai(p), A2{p) ) will be obtained by using the rotation matrix O as 




O 



A-ip) 
AM 



(A35) 



For a mode with p = but po ^ 0, the temporal component /lo(j5o? 0) is zero. The three 
spatial components are then integrated independently according to the action 



2e' 



(A36) 
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rrisea trajectories configurations 



0.02 706-5356 95 
0.03 695-6195 111 
0.04 605-5605" 97 

"configurations 1805, 1855, 1905, and 1955 were excluded from our analysis to avoid possible effects of a 
hardware failure as noted in 



TABLE I: QCD gauge configuration ensemble from [9] used for the measurements in this work. Con- 
figurations used for measurements are separated by 50 monte-carlo trajectories. Parameters for the 
ensembles are listed in Tab. |lll 



/3qcd 


0.8 


Lattice size 


V = 16^ X 32, L^=12 


domain-wall height 


M5=1.8 


Sea quark masses 


0.02, 0.03, 0.04 


l/a 


1.691(53) GeV from mp=770 MeV 
1.688(21) (IJ^^) GeV from ro=0.5 fm 



TABLE II: Simulation parameters in the generation of two-flavor dynamical QCD configurations, and 
the lattice spacing obtained in the pure QCD simulation. For details, see Ref. [8|]. 



/3em = 2/e2 


2 


e 


V^47r/137(physical value), 0.6, 1.0 


domain-wall height 


M5=1.8 


l/a 


1.659(37) GeV from mp=775.8 MeV 



TABLE III: Simulation parameters in this work. The same value of M5 as the sea quarks is used for 
the valence quarks. 
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rriu m<i = 0.015 0.02 



0.025 



0.03 



0.035 



0.04 



rUsea = 0.02 



0.015 0.2580(27 

0.02 0.2760(25 

0.025 0.2930(24 

0.03 0.3092(23 

0.035 0.3247(23 

0.04 0.3396(22 



0.015 


0.2633(22) 


0.2806(21 


0.02 


0.2808(21) 


0.2970(20 


0.025 


0.2974(20) 


0.3127(19 


0.03 


0.3132(19) 


0.3279(18 


0.035 


0.3285(19) 


0.3425(18 


0.01 


0.3131(19) 


0.3560(17 



0.015 0.2659(25 

0.02 0.2836(24 

0.025 0.3004(23 

0.03 0.3165(23 

0.035 0.3319(23 

0.04 0.3467(22 



0.2759(25 
0.2927(24 
0.3088(23 
0.3242(22 
0.3390(21 
0.3534(21 



0.2835(24 
0.3002(23 
0.3161(23 
0.3314(22 
0.3462(22 
0.3604(21 



0.2928(24 
0.3087(23 
0.3240(22 
0.3387(21 
0.3530(20 
0.3668(20 



0.2972(20 
0.3127(19 

0.3276(18 
0.3422(18 
0.3562(17 

0.3(399(17 



0.3002(23 
0.3160(23 
0.3312(22 
0.3459(21 
0.3601(21 
0.3739(21 



0.3090(23 
0.3240(22 
0.3386(21 
0.3528(20 
0.3665(20 
0.3798(19 



.03 



0.3130(19 
0.3277(18 

0.3421(18 
0.3560(17 
0.3696(17 

0.3829(10 



.04 



0.3162(23 
0.3312(22 
0.3458(21 
0.3599(21 
0.3737(21 
0.3870(20 



0.3244(23 
0.3388(21 
0.3528(20 
0.3664(20 
0.3797(19 
0.3926(19 



0.3281(19 
0.3422(18 

0.3561(17 
0.3696(17 
0.3827(16 

0.39.50(16 



0.3315(22 
0.3459(22 
0.3599(21 
0.3736(21 
0.3868(20 
0.3998(20 



0.3392(22) 
0.3530(21) 
0.3665(20) 
0.3797(19) 
0.3926(19) 
0.4052(18) 



0.3427(19) 
0.3563(17) 

0.3696(17) 
0.3827(16) 
0.3955(16) 

0.1080(15) 



0.3463(22) 
0.3601(21) 
0.3736(21) 
0.3868(20) 
0.3997(20) 
0.4122(19) 



TABLE IV: Summaxy of ud pseudoscalar 
point correlation functions. Fit range is 9 



meson masses obtained from fits to the pseudoscalar two- 
< f < 16 in each case, a = 1/137. 



rrisea ruu = 0.015 0.02 0.025 0.03 0.035 0.04 
0.02 0.2575(27) 0.2923(24) 0.3236(22) 0.3524(20) 0.3794(19) 0.4049(18) 
0.03 0.2628(22) 0.2966(20) 0.3273(18) 0.3558(17) 0.3825(16) 0.4078(15) 
0.04 0.2654(25) 0.2997(23) 0.3308(22) 0.3596(21) 0.3865(20) 0.4119(19) 

TABLE V: Summary of neutral uu pseudoscalar meson masses obtained from fits to the pseudoscalar 
two-point correlation functions. Fit range is 9 < t < 16 in each case, a = 1/137. 
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rrisea rud = 0.015 0.02 0.025 0.03 0.035 0.04 



0.02 0.2566(27) 0.2913(24) 0.3226(22) 0.3514(20) 0.3783(19) 0.4037(18) 
0.03 0.2620(22) 0.2957(20) 0.3263(18) 0.3547(17) 0.3814(16) 0.4066(15) 
0.04 0.2645(25) 0.2988(23) 0.3298(22) 0.3585(21) 0.3854(20) 0.4108(19) 

TABLE VI: Summary of neutral dd pseudoscalar meson masses obtained from fits to the pseudoscalar 
two-point correlation functions. Fit range is 9 < t < 16 in each case, a = 1/137. 



rrisea rriu = 0.015 0.02 0.025 0.03 0.035 0.04 

0.02 0.3524(22) 0.3658(20) 0.3788(20) 0.3916(19) 0.4042(18) 0.4164(18) 

0.03 0.3557(18) 0.3688(17) 0.3818(16) 0.3945(16) 0.4070(15) 0.4192(15) 

0.04 0.3594(22) 0.3728(21) 0.3859(20) 0.3987(20) 0.4112(19) 0.4235(19) 

TABLE VII: Summary of charged us pseudoscalar meson masses obtained from fits to the pseudoscalar 
two-point correlation functions, = 0.0446. Fit range is 9 < t < 16 in each case, a = 1/137. 



rrisea nid = 0.015 0.02 0.025 0.03 0.035 0.04 0.0446 

0.02 0.3512(22) 0.3645(20) 0.3775(20) 0.3903(19) 0.4028(18) 0.4150(18) 0.4261(17) 

0.03 0.3546(18) 0.3677(17) 0.3806(16) 0.3933(16) 0.4057(15) 0.4179(15) 0.4289(15) 

0.04 0.3582(22) 0.3715(21) 0.3846(20) 0.3973(20) 0.4098(19) 0.4220(19) 0.4330(19) 

TABLE VIII: Summary of neutral ds and ss pseudoscalar meson masses obtained from fits to the 
pseudoscalar two-point correlation functions, rus = 0.0446. Fit range is 9 < i < 16 in each case, 
a = 1/137. 
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meson mgea fUy^ = 0.015 0.02 0.025 0.03 0.035 0.04 

ud 0.02 0.5305(60) 0.5449(48) 0.5586(40) 0.5724(35) 0.5865(31) 0.6010(29) 

ud 0.03 0.5494(79) 0.5614(63) 0.5741(53) 0.5873(47) 0.6010(42) 0.6151(39) 

ud 0.04 0.5695(132) 0.5809(104) 0.5930(87) 0.6059(74) 0.6193(65) 0.6330(58) 

uu 0.02 0.5304(60) 0.5449(47) 0.5586(40) 0.5725(35) 0.5867(31) 0.6012(29) 

uu 0.03 0.5498(78) 0.5617(62) 0.5743(53) 0.5876(47) 0.6013(42) 0.6155(39) 

uu 0.04 0.5691(130) 0.5807(104) 0.5930(86) 0.6060(74) 0.6195(65) 0.6333(57) 

dd 0.02 0.5300(60) 0.5444(48) 0.5580(40) 0.5717(35) 0.5858(31) 0.6003(29) 

dd 0.03 0.5493(79) 0.5612(63) 0.5737(53) 0.5868(47) 0.6005(42) 0.6145(39) 

dd 0.04 0.5691(132) 0.5804(105) 0.5925(87) 0.6053(74) 0.6187(65) 0.6324(58) 

TABLE IX: Summary of degenerate vector meson masses obtained from fits to the vector two-point 

correlation functions, averaged over polarizations. Fit ranges are 5 < t < 16, 6 < i < 16, and 
7 < t < 16, for rwsea = 0.02,0.03, and 0.04, respectively, a = 1/137. 



a„ by X^/dof 


Cv dy X^/dof 


0.461(11) 2.04(17) 0.20(89) 


0.00028(20) -0.0015(31) 2.9(34) 



TABLE X: The results of the chiral extrapolation of mp± and mp± — nipO. 



fit range 6 Sq S+ 6- 5 sea ^mres X^/dof 

0.015-0.0446 4.62 (18) xlO-^ 0.0080 (12) 0.01129 (24) 0.01746(33) - 6.8 (10) xlO-^ 1.7(1.2) 

4.45 (56) xlO-^ 0.0080 (12) 0.01132 (23) 0.01741(29) 2.5(8.4) x lO"^ 6.8 (10) xlO"^ 1.7(1.2) 
0.015-0.03 4.85 (21) xlO"^ 0.0077 (20) 0.01059 (32) 0.01696(40) - 7.9 (14) xlO"^ 1.4(1.4) 

6.46 (86) xlO"-* 0.0077 (20) 0.01048 (32) 0.01701(40) -0.0028 (15) 7.9 (14) xlO'^ 0.14(25) 

TABLE XI: The results of the NLO fit to the meson mass squared differences. The first column refers 
to the quark mass range, both sea and valence, used in the fit. 



fit range Bq L5 — 


2L8 


L4. — 2^70 


xVdof 


0.015-0.0446 2.172(67) 5.7 (36) 


xlO"^ 


-0.99 (28) xlO""^ 


0.91(30) 


0.015-0.03 2.14(13) -8.6 (53) 


xlO-5 


-6.4 (77) xlO-5 


0.08(6) 



TABLE XII: The results of the a = NLO fit to the meson masses. The first column refers to the 
quark mass range, both sea and valence, used in the fit. 
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fit range m„ rud mud rUg mu/md ms/ruud 



0.015-0.0446 2.96 (13) 5.47 (11) 4.21 (12) 122.5 (27) 0.540(14) 29.08(30) 
0.015-0.03 3.02 (27) 5.49 (20) 4.25 (23) 119.5 (56) 0.550(31) 28.10(38) 

TABLE XIII: Light quark masses evaluated for physical meson masses. All values are in MeV and 
given in the MS scheme at renormalization scale /x = 2 GeV. The first column ( "fit range" ) refers to 
the range of bare quark masses used in the meson mass-squared fit. a = 1/137. 
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FIG. 1: One of the pseudoscalar mass-squared splittings computed for e = 1 (circles) and also averaged 
over e = ±1 (triangles). The latter has dramatically reduced statistical error (and is shifted for clarity). 
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FIG. 2: The quark mass dependence of charged vector meson mass and its chiral extrapolation using 
the nisea = rnval data and Eq. ([38j) . 
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FIG. 3: The quark mass dependence of the p^-p^ mass spHtting. The dashed Une is a fit to Eq. (j39p . 
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FIG. 4: Same as Fig. [3] except the minimum distance in the mass fits is t = 5 for all cases. 
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FIG. 6: The pseudoscalar mass-squared splitting for m^ea = 0.02. Each point corresponds to the 
mass-squared computed with a = 1/137 minus the same quantity computed with a = 0. The labels 
in the figure correspond to the charges of the quarks, i.e, ud means Qu = 2/3 and = 1/3. 
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FIG. 7: The pseudoscalar mass-squared splitting dependence on the sea quark mass for a representative 
case. Only splittings of degenerate mesons are shown. 
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